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, We consider the Cooper-problem on a lattice model including onsite and near-neighbor interac- 

' tions. Expanding the interaction in basis functions for the irreducible representation for the point 

Q\ \ group Civ yields a classification of the symmetry of the Cooper-pair wave function, which we calcu- 

late in real-space. A change of symmetry upon doping, from s-wave at low filling fractions, to d x 2_ y 2 
at higher filling fractions, is found. Fermi-surface details are thus important for the symmetry of the 
superconducting wave function. Symmetry forbids mixing of s-wave and d-wave symmetry in the 
Cooper-pair wavefunction on a square lattice, unless accidental degeneracies occur. This conclusion 
04 ' also holds for the selfconsistent treatment of the many-body problem, at the critical temperature T c . 

CI Below T c , we find temperatures which are not critical points, where new superconducting channels 

open up in the order parameter due to bifurcations in the solutions of the nonlinear gap-equation. We 
(— i | calculate the free energy, entropy, coherence length, critical magnetic fields, and Ginzburg-Landau 

parameter k. The model is of the extreme type-II variety. At the temperatures where subdomi- 
nant channels condense, we find cusps in the internal energy and entropy, as well as as BCS-like 
discontinuities in the specific heat. The specific heat anomalies are however weaker than at the true 
fS |i superconducting critical point, and argued to be of a different nature. 
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I. INTRODUCTION 



The non-perturbative effect on the ground state wave-function of an electron gas with arbitrarily weak effective 
attraction between the quasi-particles on the Fermi surface was first demonstrated by Cooper in his essentially exact 
solution of the corresponding two-body problem Regardless of the strength of the interaction, two electrons 
interacting attractively on opposite sides of the Fermi-surface in an otherwise inert Fermi-gas form a bound state, 
the Cooper-pair. This simple calculation alone suffices to yield precisely the correct non-analytic dependence of the 
binding energy of the Cooper-pair on the dimensionless coupling constant A later found by solving the full selfconsistent 
problem Q). Cooper solved the problem for the simple "jellium" model of a metal, where the Fermi-sea was taken to 
be spherical, and the interaction between the two extra electrons of s-wave symmetry. Of primary interest was the 
two-body spectrum, the Cooper-pair wave-function in k-space under such circumstances being a trivial constant. 
' In this paper, we reconsider this simple problem on a square lattice. The tight-binding band structure includes 
, nearest and next-nearest neighbor hopping, while the two-body term in the Hamiltonian includes an onsite repulsive 
f^i • Hubbard-term, a nearest neighbor effective electrostatic interaction, which may be taken to be attractive, and also a 
t^- \ next-nearest neighbor electrostatic interaction. The problem now includes two additional non-trivial features: i) The 
ON i Fermi-surface is no longer spherical and one has to work directly in k-space rather than transforming the problem 
to energy-space, ii) The interaction between the quasiparticles no longer has simple s-wave symmetry, but rather 
may be expanded as a bilinear combination of basis functions for the irreducible representation of the point group 
of the 2D square lattice, Ci v . This leads to the possibility of a number of interesting effects. In addition, we use a 
selfconsistent scheme to calculate the superconducting gap, thermodynamic quantities, and temperature dependence 
■ of critical magnetic fields of this phenomenological lattice fermion model. 

This paper is organized as follows. In Section II, we define the model to be considered. The method of calculating 
the Cooper-pair wavefunction is presented in Section III, while specific numerical results pertaining to this quantity 
are presented in Section IV. More detailed analytical results on the binding energy of the Cooper-pair are given in 
Section V. In Section VI, we present the results for thermodynamic quantities and critical magnetic fields from a self- 
consistent scheme for a gap-function with several symmetry-channels, belonging to various s-wave and d-wave channels. 
In Section VII, we give a discussion of the specific heat anomalies one may expect in such a model. We emphasize that 
throughout this paper, the superconducting order parameter is not a vector order parameter, but assumed to be a spin- 
singlet scalar complex order parameter such as is believed to describe conventional low-temperature superconductors 
and high-T c cuprates. 



II. THE MODEL 



The model we consider is an extended Hubbard-model on a square lattice defined by the Hamiltonian 
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Here (i, j) and denote nearest neighbor and next-nearest neighbor couplings, respectively, t, and t' are corre- 

sponding hopping matrix elements, and \x is the chemical potential. U is an onsite repulsion term, while V and W 
are effective electrostatic Coulomb matrix-elements between nearest and next-nearest neighbors, respectively. When 
viewed as an effective interaction term within a one-band model, as reduced from a multiband-band model, these 
terms may be attractive &-^]- In this paper we simply take them as effective attractions without further justification. 

After introducing a plane-wave basis and performing a standard BCS truncation of the interaction piece of the 
Hamiltonian , it takes the usual form 

H = y^e r d cr +V y rp c+ c+ c p , cd t , (2) 

Z—i k k,a k,a k,k' j^f -k,l ~k ' ,1 fe',T' v ; 

k,(T k,k' 

where we have denned, after an appropriate redefinition of the zero-point of energy and a rescaling of t' and fi 

Ej: = — 2t[cos(fe x ) + cos(ky) — It! cos(k x ) cos(k y ) — (2 — 2i' — n)] . (3) 

The choice of such a quasiparticle dispersion is obviously motivated by its relevance as a simple means of modeling the 
quasi-particle band crossing the Fermi-surface of the high-T c cuprates ||. Note that in such a context, the inclusion 
of the i'-term is crucial; a bipartite lattice with nearest-neighbor hopping only, is inconsistent with the observed 
Fermi-surfaces in the high-T c cuprates, where \t'\ w |i|/2 ||. (This also has motivated the choice t' = 0.45t in our 
numerical calculations). The importance of including this term in correctly interpreting experiments, has recently 
been strongly emphasized 

In the above truncation of the interaction term, the interaction is assumed to be operative between fermion spin- 
singlets on opposite sides of the Fermi-surface, the inert and rigid Fermi-sea merely giving rise to Pauli-blocking 
factors. With the interactions given in Eq. [I], it is readily shown that p is given by 
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V U,P = E A " B v$)> ( 4 ) 

where Ai = U/2, and A2 = A4 = V, A3 = A 5 = W. We have also found it convenient to introduce the simple, but 
sufficient, subset of basis functions {B n (k)} for irreducible representations of the symmetry group C± v of the square 
lattice, Bi(k) = -^=, B 2 (k) = ^[cos(fc x ) +cos{k y )}, B 3 (k) = -^[cos(fc K ) cos(fcy)], B 4 (k) = -^=[cos(k x ) ~ cos(k y )], 

Bs(k) = ^^[sin(fc x ) sin(fcj,)], where N is the number of lattice sites. Inclusion of longer ranged interactions will in 
general require an augmentation of this subset, but any finite ranged interaction will yield a separable potential. 



III. THE COOPER PAIR WAVE FUNCTION 



We define a two-particle state for the non-interacting case, i.e. U = V = W = 0, by \k, a; — fc, — c)o obeying the 
Schrodinger equation 

H \k, a; -k, -a) = 2£^|fc, a; -k, -a ) , (5) 

where H denotes the Hamiltonian of the free particles, and eg is given in Eq. ||. Note that in this notation, the 
hopping between next nearest neighbors is given by the matrix element, 4t * t', and we limit ourselves to situations 
where 2\t'\ < I, such that the bottom of the band is located at the Brillouin-zone center. 

The problem we will consider is a simplification of the one posed by H = Hq + H-^. We imagine that we have a rigid 
Fermi-sea of non-interacting electrons with a spectrum given by Eq. |^. To this inert Fermi-sea we add two electrons 
on opposite sides of the Fermi-sea which interact with the matrix element Vk,w ■ This interaction term scatters a pair 
of electrons in the state \k, <r; —k, — er)o to the state \k' , a; —k', —<t)q. The exact two-particle state for the two extra 
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electrons, for which k is no longer a good quantum number, is defined by |1, 2). This state is expanded in two-particle 
plane-wave states as follows 

|1,2)= J2 H,„ \k,<r;-k,-<T) , (6) 

k>kF ,cr 

such that the problem is reduced to one of determining the Fourier-coefficients . The exact two-particle state for 
this problem obeys the Schrodinger equation: 

(H + H int )\l,2) =E\l,2). 

Upon inserting the expansion in plane-wave states, and projecting onto plane-wave states, we obtain upon using Eq. 
^, an integral equation for the expansion coefficients 

£ B v (k') X v H B S) = (E- 2s g) op, 

r\ k>kp 



which immediately yields an expression for the expansion coefficients 

a k> = E _ 2 > £k>e kF . (7) 

Here B r) (k) is one of the five basis functions required to expand the interaction. Finding the wave function, thus 

amounts to finding the eigenvalue E and the amplitudes A n . Inserting the expression obtained above for and V?z, 
into the Schrodinger equation, we obtain coupled algebraic equations for the amplitudes, A v as follows 

Ev^(^')E^ E Bv eI b 2 ^ ] = Y, A v> B Ak')- (8) 

77' rj k>k F k 7}' 



The functions B r) (k) are linearly independent, and a comparison of coefficients therefore yields the following set of 
equations for the amplitudes A v 

J2 A i V D m' = A V- ( 9 ) 

n 

The result can be written as a vector equation, (T — I) A = 0, where T ri ' n — X^D^i and I denotes the identity 
matrix. A nontrivial solution exists if and only if the system determinant vanishes |T — I| = 0, which determines the 
eigenvalue E, and thus in turn the eigenvectors. 

Note that this analysis shows that at this level, in general one cannot obtain a Cooper-pair wave function with a 
mixed (s, G?)-symmetry. This is a consequence of the fact that the matrix -D w ' is block-diagonal in the s-wave and 
d-wave sectors due to the fact that the "pair-susceptibility" Xfc = 1/(E — 2eg) transforms as a function with s-wave 
symmetry expandable in the functions B\, Bi and B3. In Section VI, we show that this conclusion holds for the 
mean- field gap-equation that results from a self-consistent solution to the full problem at the critical point T c . The 
wave-function may thus have one of the following forms 



^2 A v B v( k ) Xfe ; £ fe > £ k F , 

ri=l 

5 

g = E A *i B i$) X£ ; £k > e kF . (10) 



»j=i 
5 
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"k 

■q=4 



All the amplitudes are zero for < Sk F - Since the D-matrix is block diagonal, we may find two different 
eigenvalues, one for each irreducible representation. The correct eigenvector corresponds to the lowest eigenvalue E. 
Furthermore, an immediate consequence of the fact that = ^ k>kF A I( a k B v (k) is that A n 7^ if and only if 
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A ?) 7^ 0. Specifying the X v 's thus immediately determines which of the basis functions B v (k) that may contribute to 
the Cooper pair wave function. 

The Cooper pair wave function in real-space is calculated by applying the inverse lattice Fourier-transform to , 
defined by 

nk— Ink— 1 

= £ ¥ i|((fe ' fe>(i ' j>) ' (id 

fc 1= k 2 =0 

where f — is the relative distance between the electrons comprising the Cooper-pair, k — (fei,^). 



IV. NUMERICAL RESULTS 



In Fig. |j] (a)-(f) we have plotted \i^(i,j)\ 2 of the Cooper pair wave function given by Eq. ( pi] ) in the case of no 
onsite interaction, U — 0, and nearest neighbor attraction, V = — 0.75*i. Throughout this discussion, we set W — 0. 
The hopping matrix element t = 0.10 eV. At low filling fractions, the wave function displays s-wave symmetry which 
can be recognized by the peak present at the point (0,0). Electrons situated at the same site forming a pair cannot 
be in a d-wave state. An example of a d-wave state is shown in Fig. |l|(c). This state is in fact found when increasing 
the filling fraction to n = 0.17, for the same value of the Coulomb-parameters. The pairing state has changed its 
transformation properties due to doping. Further increasing the filling does not change the symmetry, but Figs. 
|l|(c)-(f) show that the pairing state broadens as doping is increased. The region where \4>{i,j)\ 2 ^ can roughly be 
interpreted as the size of the Cooper pair, i.e. the coherence length, which is thus seen to increase with increasing 
filling fraction. This is verified directly by calculating the quantity J2i j I* — j\ 2 l^ih j)\ 2 / J2i j IVK^i j)\ 2 ■ 

To further investigate at which doping level the symmetry change occurs, we have plotted the eigenvalue E of the 
two-particle Schrodinger equation, in Fig. |^ as a function of doping. The eigenvalue E is determined by numerically 
solving |T — I| = 0. The Fig. shows that for small filling fractions, i.e. n < 0.1, the eigenvalue for s-wave pairing is 
lower than the corresponding value for d-wave pairing, and hence s-wave pairing is therefore energetically favorable. 
At n « 0.10, we see that the two energy curves intersect, and for n > 0.10, d-wave pairing is thus expected to be 
energetically favorable. 

We have done similar eigenvalue calculations as above in the case of an onsite repulsion, U = 1.0 * i. The results 
are shown in Fig. ^. Since an onsite interaction does not affect d-wave pairing, the curve that shows the d-wave 
eigenvalues is the same as shown in the scenario of no onsite interaction. The s-wave curve has however changed, and 
appears to be shifted to higher energies, which means that the change in the symmetry of the wave-function occurs at 
a lower filling fraction, in this case n ~ 0.08. The results presented in Fig. indicate that increasing onsite repulsion 
favors d-wave pairing, as one would expect. The extended s-wave component B2 of the Cooper-pair wave function 
will inevitably have a finite on-site component, although it avoids the hard core to a considerable extent. 

Figs, ^(a)-(f) show the results for \ip(i,j)\ 2 as a function of doping for U = 4.0 * t. This onsite repulsion is 
large enough to suppress pairing in any of the s-wave channels for all fillings that we have considered, in the range 
n G [0.06, 0.85] . As filling is increased, one starts seeing a considerable broadening of the wavefunction with increasing 
filling fraction. 

Figs. H (a)-(f) show the Cooper-pair wave function with n — 0.06 at six different values for the onsite repulsion 
U. In Fig. |^(a), the onsite repulsion {7 = 0, and the s-wave solution is the energetically favorable (see Fig. ||). In 
Fig. |](b) a weak onsite repulsion U = 0.5 * t is applied. Note that for this case, we obtain the somewhat counter- 
intuitive result that, since j)| 2 ; \i — j\ =0 evidently has increased upon increasing the onsite repulsion, the weak 
increased repulsion effectively promotes attraction between electrons. Similar effects have been seen in exact many- 
body calculations on strongly correlated ID lattice fermion models pj . We will comment on this result later. Figs. 
|(c)-(f) display the wave function when the onsite repulsion is increased in the range U/t — (0.7, 1.0, 1.30, 1.45) and 
it is evident that between U = 1.30 * t and U = 1.45 * t, the wave function has changed its transformation properties 
to d-wave in order to completely eliminate the effect of the hard core. Increasing the onsite repulsion further will not 
affect the d-wave function, as already mentioned. 

To show in more detail for which values of U the wave function changes its transformation properties, we have 
calculated the energy eigenvalues as a function of onsite repulsion for a filling fractions n = 0.06. The results arc 
shown in Fig. |6[ The trends in the results clearly show that in lattice models with four-fold symmetry, s-wave 
superconducting pairing is a low-filling effect, while d-wave pairing practically always wins in situations close to 
half-filling. 
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V. ANALYTICAL RESULTS 



It is instructive also to perform a simplified analytical treatment of the Cooper-problem. A major simplification 
results by rewriting D vv > (E) as follows 



D m , (E) =1 de^ S(e - e s ) 9(e s ) B ^_ B ^ = J de^M 9(e) 



k 



N W I de ^— ^ - N W F(—), (12) 



where we have introduced the binding-energy A = 2fi~ E, lo c is an upper band-cutoff and the function F is defined as 
F(x) = \ ln(l + 1). In addition, we have introduced the projected densities of state N, pl /(e) = J2k W 
e^). The constants N vri i in Eq. [l^ are chosen appropriately to get the correct value for the /c-space sum. With our 
choice of coupling constants, the matrix block-diagonalizes into a 3 x 3 matrix in the s-wave sector, and a 2 x 2 
matrix in the d-wave sector. For simplicity, in the following we will set W = 0. It will also be convenient to introduce 
the following definitions 

v = ^ (N n N 22 - N 12 2 ); 7=1 (|jV u + VN 22 ). (13) 

With these definitions, we get the following secular equation for F (and hence A) in the s-wave sector v F 2 + 2j F+l = 
0. This is a necessary, but not sufficient, condition for finding a non-zero Cooper-pair wave function in the s-wave 
sector. Note that we must in principle distinguish between four cases: i) v < 0, 7 > 0, ii) v < 0, 7 < 0, iii) 
v > 0, 7 > 0, iv) v > 0, 7 < 0. Here, we will discuss the case v < 0, which should be realized for appropriate 
band- fillings. Furthermore, it is clear that in the large-?/ limit, 7 > 0, while 7 < when U = 0. Hence, we have 



F _ l7l + y/|7l 2 + H 2 J_ TJ ^ \V\ N22 
M 2\' U> 2 



p H7l+ Vl7| 2 + H 2 1 TT , \V\ N 22 

which defines the dimcnsionlcss coupling constant X s . It follows that when U = 0, we get F = 2/\V\ N 22 while for 
U = 00, we find F = 2/(\V\ {N 22 — N 21 /Nu). Note that, within the s-wave sector, even an infinite hard-core repulsion 
does not suffice to destroy superconductivity. This is mainly due to the fact that the extended s-wave piece of the 
wave-function avoids the hard core. What will ultimately determine the symmetry of the superconducting wave- 
function is therefore the competition with the d-wave sector. In the d-wave sector, we obtain the secular equation 
1 - 2X d F = 0. 

In either case, to have a bound state form, it is required that F is real and positive, equivalently the effective 
dimensionless coupling constants A s and Xd must be positive. The only acceptable solution to the secular equation in 
the s-wave sector is F = 7^-; A s — — ttttH — r, while for the d-wave sector we obtain F — ttI— ; A^ = — , 2 f. c , . . 

2A S ' 3 cxp(l/A s ) — 1' 2Arf ' u cxp(l/Ad) — 1 

The binding energy, in the s-wave and d-wave sectors, are determined by the effective coupling constants A s and Xd, 
respectively. The system chooses the symmetry of the superconducting wave-function which gives the largest binding 
energy, equivalently the largest A. 

To gain some insight into what filling fractions one might expect one or the other coupling constant to dominate, it 
is instructive to exhibit N 22 (e) and 7V44 (e) as a function of the energy. They are shown in Fig. 0, along with the 
single-particle density of states N(e) — J^g 5(e — et), which we have been able to evaluate analytically. This is useful 
since it provides intuition about the singularities that also appear in the projected densities of states, and which we 
have to consider when the stability of the various superconducting channels arc investigated. Using Eq. ^|, we find 
straightforwardly that the single-particle density of states is given by [12] 



m = 1 » k(,HM5); | ((I + „)/2| < 1 

2 tx z \t\ \J\ — pa VU 1— pa / 



1 1 . /8 VI -pa 



ln ,«VWaV {{a + m<1 
V a — o / 



2 7T 2 \t\ VI ~ 9 a ""V \a-p\ 
= 0; |(o + p)/2| > 1, (15) 
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where we have defined p — —21/, a = e/2t, and K(fc) is a complete elliptic integral of the first kind, with modulus k 
fl3|| . The result is valid for \p\ < 1, the requirement for having a band-minimum at the center of Brillouin-zone. (It is 



quite remarkable that the approximate logarithmic expression in Eq. 15 reproduces the exact expression excellently 
over the entire band. The reason is that the main effect of the elliptic integral K is to produce a narrow logarithmic 
singularity at a = p, while the remaining variation in N(e) is due to the prefactor.) The density of states has a 
logarithmic singularity N(e) = (l/2n 2 \t\y/l — p 2 ) [ln(8/|a — p\) + ln(^/l — p 2 )] at e — 2 t p, a finite cusp at the 
lower band-edge, N(e = 2t(—2 — p)) = (l/47r|i|)/(l + p), while the density of states close to the upper band-edge 
has a smooth behavior with a value at the upper band-edge given by N(e = 2t(2 — p)) = (l/47r|i|)/(l — p). Were we 
to reverse the sign of t' , we would get a cusp at the upper band-edge, while the density of states at the lower edge 
would be smooth. The singularities that appear in the projected densities of states must come from the singularities 
in N(s). We have not been able to reduce the integrals for the projected densities of states to useful expressions in 
terms of elliptic integrals, and have simply calculated them numerically using the tetrahedron algorithm p"4[ . 
The peaks in N22 (e) originate from energies close to the bottom of the band, essentially identified with the cusps in 
N(e), whereas the peaks in the A^fe) come from energies identified with the logarithmic singularities in N(e). As 
soon as the filling fraction increases from the bottom of the band, the major singularities in iV22(e) are eliminated from 
the integrated projected density of states N22, and the contributions to JV44 coming from the logarithmic singularities 
in N4i(e) at e = 2tp will dominate. Hence, Ad and thus d-wave pairing will always dominate close to half- filling, 
while s-wave pairing will win out close to the bottom of the band. Note for instance in the case where t' = 0, that 
^22(2) = in the middle of the band, while A^fe) has a weak logarithmic singularity. When t' the situation 
changes somewhat, and both projected densities of state are finite at half filling. Nonetheless, as long as t' < 0.5, such 
that the bottom of the band is located at the zone-center, N44 > N22 for filling fractions such that the Fermi-surface 
is close to the middle fo the band. For non-pathological situations, d-wave pairing will in general dominate close to 
half-filling. The conclusion holds for any pairing kernel with a predominant .E^-symmetry. 

A somewhat counter-intuitive result obtains if we introduce a weak repulsive coupling U in the problem. In our 
notation, this means that in Eq. [l4|, we may expand F in Eq. [l4]in powers of M/Mi m order to find the corrections 
to A s to leading order in U/V. The leading order correction in U is given by 

^~h\ = ^- r P\. (16) 

Thus, we see that an enhancement of X Sl and hence an enhancement of superconductivity, comes about initially as U 
is switched on, irrespective of the sign of U. On the other hand, if we consider A s at U — 00, it is clearly less than X s 
at U = 0, 

X S (U = 0) - X S (U = 00) = M ^! > 0. (17) 

l iVn 

It follows that there must exist an intermediate value of U which maximizes the binding energy in the s-wave sector. 
Since Xd is independent of U altogether, we obtain the counter-intuitive result that there exists a regime of filling- 
fractions, where one may obtain a "switch" of the symmetry of the superconducting wave function from d-wave to 
extended s-wave upon increasing a weak onsite repulsion. This has much in common with the result we observed in 
the numerics, that an increase in U for n — 0.06 produced a more tightly bound Cooper-pair. While the analytical 
results quoted above are approximate, the numerical results are quite remarkable in that they, within the Cooper- 
scheme, entail no further approximations. We also note that similar results have been found for the zero-temperature 
superconducting gap in self-consistent mean-field calculations Jig ], We now comment on why we may obtain such a 
counter- intuitive result in our numerical calculations. 

A finite interaction V < implies that when U = 0, there certainly will exist a solution to the Schrodinger equation 
such that A2 7^ 0, i.e. we are guaranteed pairing in the extended s-wave channel for low enough band-fillings, as 
argued above based on our intuition on the projected densities of states. Increasing U slightly from zero, implies that 
one in principle is coupling in another component of the wave-function, A±, since Ai ~ Ai ~ U. However, A\ and 
A2 may enter in the wavefunction with relative phases of it. The coupling between Ai and A2 in the gap-equation 
may take advantage of this and essentially produce an "attraction from a repulsion" by appropriately twisting the 
relative phases of the amplitudes entering in the order parameter. Note that such a mechanism cannot work with 
an order parameter containing only one symmetry-channel component, since the overall phase of the wavefunction 
is irrelevant for the spectrum. In the present situation with two s-wave components, however, the extended s-wave 
amplitude is able, for small values of U, to "drag" the isotropic s-wave component of the order parameter along, as 
a consequence of the "conversion of repulsion to attraction" -effect. Continuing to increase the onsite repulsion to 
larger values U/t » 1, energetics dictates that the on-site component of the order parameter eventually must vanish, 
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thereby reducing the indirect attractive effect of the onsite-repulsion. As a consequence, the pair-wave function will 
again spread out in real space, as observed in our numerics. 

Note that the above results are not sensitive to the basic pairing mechanism operative in the underlying microscopic 
theory. All that enters at this level, are competing channels in the pairing kernel with isotropic and extended s-wave 
symmetries, as well as d-wave symmetries. Hence, for instance anti- ferromagnetic spin-fluctuations who some 
believe to be relevant in producing superconductivity in the high-T c cuprate CitC^-planes, would fit into the scenario 
Eq. H for the pairing kernel. Nonetheless, the above "attraction from repulsion" effect within the isotropic and 
extended s-wave sector is not particularly relevant to the physics of the high-T c oxides, which show superconductivity 
only in the vicinity of half-filled bands, where we find that d x 2_ y 2-w&ve pairing almost invariably wins out. 



VI. SELFCONSISTENT ANALYSIS 

In this Section, we consider the full nonlinear mean-field gap-equation. Starting from the full Hamiltonian, Eq. [I], 
one performs a standard BCS-truncation of the interaction term, and a further anomalous mean-field decomposition 
to obtain the gap-equation 

k< 

V £,£>=E \B n {k)B n {k') 
n 

XT = tanhf — k 



^ = V £ £ + |A * 12 - (18) 

Here, Sj: is the normal state dispersion relation, and f3 = 1/fceT. Given the form of Vgc, , it is clear that the 
gap-function must be expandable in the basis functions for the irreducible representations of C<± v in the following way 

Ag = J2 A n B v(k), (19) 
v 

where A^ are amplitudes that must be found selfconsistently. Inserting such an Ansatz back into the gap-equation 
and equating coefficients of the linearly independent functions {_B^(fc)}, one finds the following coupled, nonlinear 
algebraic equations for the gap-amplitudes A,, 

A„ = J2 A v' M vn' 
n' 

M w = -X v ^ B v (k) B v ,(k) XS . (20) 

k 

Notice that the pair-susceptibility xj: transforms as an s-wave function. Hence, the matrix M vri i must be block- 
diagonal by symmetry, in the s-wave and d-wave sectors. Including the five basis functions we have used in this paper, 
the matrix thus block-diagonalizes into a 3 x 3 and a 2 x 2 matrix. 

Consider first the case T — T c , where the equations are linear. The various sectors must be solved separately for 
the s-wave gap and the d-wave gap, respectively. The physically relevant solution is the one with the largest value of 
T c , while the other solution is unphysical. The alternative is a scenario with two distinct superconducting transitions 
in zero magnetic field, which is not acceptable. Therefore one can never get a mixing of s-wave and d-wave symmetry 
on a perfectly square lattice at the critical temperature. The possible exception is the one where coupling constants 
in combination with filling fractions have been fine-tuned in such a way that the two sectors yield the same T c . 

Below T c , the situation is changed. The matrix M nv / is still block diagonal, but there is nonetheless coupling between 

the s-wave and d-wave sectors, since the pair-susceptibility depends on the entire gap A^ = A J; B n (k). As a 
consequence, at T = T c the dominant channel becomes superconducting, while more components are added to the gap, 
at distinct temperatures, as the temperature is lowered. This means that one is getting superconducting condensation 
in progressively more channels. Nonetheless, these subsequent additions of channels to the order parameter do not 
constitute separate superconducting critical points. They simply add more amplitude to |A^|. In mathematical terms, 
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this is equivalent to a bifurcation in the solution of the nonlinear gap-equation. Such bifurcations have recently been 
studied for the BCS-gap equation [^7| . 

We have solved the coupled mean-field gap equations numerically, to obtain the superconducting gap, and hence 
various thermodynamic quantities. We omit details of our calculations of the gap-function itself, suffice it to say 
that we have reproduced in detail the results of Ref. Jl5| . The free energy, internal energy, entropy, specific heat, and 
critical magnetic fields are found by using the standard expressions for the free energy F and the entropy S 

F = T,( £ k~ Xi I A £ | 2 ) - i M2 + expGMfc) + eM-PEj;)] 

k k 

S = -2k B ]T ln(l-/g) + 4 In/g], (21) 

k 

with F = U — TS, and C = TdS/dT, where U is the internal energy, and fz= + exp(/3E^)] is the momentum 
distribution function for the particle-like elementary excitations of the superconductor. 

The lower and upper critical fields B r \ and B C 2, respectively, are found from the standard expressions F(T = 
T c ) — F(T < T c ) — and B C 2 = <I>o/27r£ 2 , where we have calculated the superconducting coherence length 

from £ = vp/A. Here vp is the Fermi- velocity in the normal metallic state, while A is taken to be the absolute value of 
the gap Ag averaged over the Fermi-surface. As a sideproduct, wc obtain the Ginzburg-Landau ratio k = B c i j \f2B c \ . 

Our results exhibiting the quantities described above, are shown in Figs. [8|-|To|. As the temperature is lowered, one 
will observe cusps in |A^|, and hence the internal energy and entropy, at the temperatures where new channels are 
coupled into the gap. The specific heat will show BCS-discontinuities at all the temperatures where the new channels 
condense, in addition to the BCS-discontinuity at the superconducting transition p7| . However, the amplitude- 
fluctuations of the orderparameter are massive at these lower temperatures, and the superconducting correlation 
length is finite. This is a generic feature of superconductivity in systems with competing pairing channels, irrespective 
of whether the dominant pairing is d-wave or s-wave. 



VII. SPECIFIC HEAT 



For the case of extreme type-II superconductors, critical fluctuations will surely modify the results for the specific 
heat close to the critical temperature. The anomaly in the specific heat at the critical temperature is therefore not 
expected to be of the BCS-type at all, due to phase-fluctuations of the order parameter. As recently emphasized 
|18fl , in extreme type-II superconductors with a low superfluid density, the dominant contributions to the specific heat 
anomaly are expected to be phase-fluctuations in the order parameter, not amplitude-fluctuations as in BCS. This is 
easily seen by noting that the superfluid stiffness is related to the free energy F and a phase-twist 8<p in the order 
parameter across the system via the Fisher-Barber- Jasnow relation [fl9| 



d 2 F 
d{8<t>f 



(22) 

i=0 



On the other hand, superconductors arising from poor conductors such as the high-T c compounds, do in fact have 
a low superfluid stiffness |Ts|| . This implies that phase-fluctuations in the superconducting orderparameter are soft, 
dominating the fluctuation spectrum, while amplitude fluctuations may be neglected. One may think of the super- 
conductivity as arising out of a quantum fluid of preformed pairs. Hence, for weakly coupled layers such as those 
we have considered in this paper, the critical point is in the universality class of the 31? XY-model, and the specific 
heat should have an analagous critical anomaly, which turns out to be an asymmetric logarithmic singularity with a 
specific heat exponent a — —0.007(6) rather than the finite discontinuity of the BCS-type |po| , pi[ . This has recently 
also been observed in experiments on the extreme type-II Y BaiCuzOi cuprate with a Ginzburg-Landau parameter 
k ~ 50 |2^| . (The rather counter-intuitive result that fluctuations may convert and sharpen a mean-field like finite 
discontinuity in the specific heat into a (logarithmic) singularity at the transition, rather than smoothing over the 
asperities in the specific heat, is known from Onsager's famous solution of the 2D Ising model p3fl .) 

The weak anomalies at the lower temperatures should however be well captured by the mean-field theory, provided 
they are located outside the critical regime, which is likely to be the case. This is therefore a way of distinguishing 
intrinsic anomalies in the specific heat from those anomalies one will observe in multi-phase compounds with a 
distribution of critical temperatures. In the former case, one should observe one 3D XY-like anomaly at the critical 
temperature, and additional weaker mean-field BCS-like anomalies at lower temperatures. In the latter case, one 
should observe a number of anomalies roughly of the same order of magnitude, all of the 3D XY-type. 
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A large onsite repulsion U ultimately suppresses pairing in the isotropic and extended s-wave channels. This may 
explain why two anomalies in the specific heat so far have not been observed in zero magnetic field in Y Ba-iCu^O-j 
p2| . However one should note that an inclusion of the M^-term in Eq. [l] will give pairing in another d-wave channel, 
d xy , which could compete with d x 2_ y 2 in the presence of a very large Hubbard-[/, in the relevant doping regime, and 
hence also give two anomalies in the specific heat even in large-t/ compounds such as the high-T c cuprates. Recent 
experiments show intriguing features in microwave conductivity as well as London penetration depth well inside the 
superconducting phase, which may be consistent with the above picture [ p4| . 
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LIST OF FIGURES 

FIG. 1. The probability density ui(i,j) = \ip(i,j)\ 2 of the Cooper-pair wave function with U — at different filling fractions. 
In Figs, (a) and (b) the wave function exhibits s-wave symmetry for n = 0.06 and n — 0.12, respectively. Figs, (c) - (f) all 
exhibit d-wave symmetry of the wave function where n = 0.14 in (c), n — 0.17 in (d), n = 0.51 in (e), and n = 0.85 in (f). 

FIG. 2. The binding energy of the Cooper-pair as a function of doping, at U = 0. For n < 0.11 the binding energy of 
the s-wave pairing is energetically favorable compared to d-wave pairing. As doping increases, the d-wave pairing becomes 
favorable. 

FIG. 3. The binding energy of the Cooper-pair as a function of doping, at U = 1.0 * t. For n < 0.11 the binding energy of 
the s-wave pairing is energetically favorable to d-wave pairing. As doping increases, the d-wave pairing becomes favorable. 

FIG. 4. The probability density u)(i,j) = \ip(i,j)\ 2 of the Cooper-pair wave function with U = 4.0 * t at different filling 
fractions. Figs, (a) - (f) all exhibit d-wave symmetry of the wave function where n = 0.06 in (a), n = 0.12 in (b), n = 0.14 in 
(c), n = 0.17 in (d), n = 0.51 in (e), and n = 0.85 in (f). 

FIG. 5. The probability density co(i,j) — \tp(i,j)\ 2 of the Cooper-pair wave function with n — 0.06 at different onsite 
repulsions. For U — 0, s-wave pairing is favorable (a). For the values U = 0.50 * t,U = 0.70 * t,U = 1.0 * t and U = 1.30 * t as 
shown in (b), (c), (d) and (e), respectively, the Figs, exhibit s-wave symmetry. For U = 1.45 * t the Cooper pair wave function 
tranforms as a d-wave. 

FIG. 6. The binding energy of the Cooper-pair with n = 0.06 as a function of onsite repulsion. Since the onsite potential 
couples to an s-wave, the d-wave binding energy is not affected by increasing U. For U < 1.4 *t, s-wave paring is favored. For 
larger U's, d-wave pairing is favored. 



FIG. 7. The projected densities of states N22 and -/V44 along with the single particle density of states N(e) = S(e — eg). 

FIG. 8. Thermodynamic quantities in superconducting and normal states for n = 0.25 and U = 0. Figs, (a)- (d) show the 
entropy, free energy, internal energy and specific heat, respectively. They all exhibit a critical behavior at T ~ 90 K. In (d) 
it can easily be seen that another symmetry channel switches on at T w 60 K, giving rise to cusps in the internal energy and 
entropy, and a discontinuity in the specific heat. Figs, (e) and (f) show the coherence length and the critical fields, respectively. 
Note that the coherence length has singular behavior only at T w 90A. 
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FIG. 9. Thermodynamic quantities in superconducting and normal states for n = 0.60 and U = 4.0 * t. Figs, (a)- (d) show 
the entropy, free energy, internal energy and specific heat, respectively. They all exhibit a critical temperature at T « 90 K. 
Figs, (e) and (f) show the coherence length and the critical fields, respectively. In this case, there are no additional features 
below T c , due to the complete suppression of the instability in the competing s-wave channel as a result of the large onsite U. 



FIG. 10. The Ginzburg-Landau parameter, k, as a function of temperature for the cases: (a) n = 0.25 and U = and (b) 
n = 0.60 and U = 4.0 *t. The result shows that the superconductor model is of the extreme type-II variety. Hence, one expects 
the dominant critical fluctuations in the order parameter to be phase-fluctuations, not amplitude fluctuations. 
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